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Abstract: We perform an exploratory lattice QCD simulation of DD scattering, aimed at 
determining the masses as well as the decay widths of charmonium resonances above open 
charm threshold. Neglecting coupling to other channels, the resulting phase shift for DD 
scattering in p-wave yields the well-known vector resonance "0(3770). For m 7 r = 156 MeV, 
the extracted resonance mass and the decay width agree with experiment within large 
statistical uncertainty. The scalar charmonium resonances present a puzzle, since only the 
ground state Xco{^P) is well understood, while there is no commonly accepted candidate 
for its first excitation. We simulate DD scattering in s-wave in order to shed light on 
this puzzle. The resulting phase shift supports the existence of a yet-unobserved narrow 
resonance with a mass slightly below 4 GeV. A scenario with this narrow resonance and 
a pole at Xco(l-P) agrees with the energy-dependence of our phase shift. Further lattice 
QCD simulations and experimental efforts are needed to resolve the puzzle of the excited 
scalar charmonia. 
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1 Introduction 

Charmonia cc well below open-charm threshold DD are among the best understood hadrons. 
Their spectra and selected transition matrix elements are successfully described by lattice 
QCD simulations and QCD motivated models. Recent lattice calculations have performed 
the necessary extrapolations and considered spectra [1, 2] as well as certain radiative tran¬ 
sitions [3, 4]. For states well-below open charm threshold, the main remaining uncertainty 
is the neglect of charm-annihilation Wick contractions in lattice simulations. 

The most interesting charmonium and charmonium-like states lie near or above open 
charm thresholds. During the past decade a plethora of states that can likely not be inter¬ 
preted as conventional cc have been discovered in experiment (for a review see for example 
[5, 6]). These states have been treated theoretically making simplifying assumptions and 
reliable quantitative results for those hadrons are not available. In particular, all of the 
lattice simulations so far have ignored the strong decay of the charmonium resonances to 
a pair charmed mesons cc —)■ D^*^D^*\ which typically represents the main decay mode. 
Except in a few simulations [7-10], the effect of the threshold on the near-threshold states 
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has been neglected. The most extensive spectrum of charmonia has been obtained in sim¬ 
ulations with Nf = 2 + 1 dynamical flavors at m-,^ ~ 400 MeV [11], but the determination 
neglects the unstable nature of the states and relies on extracting the energy levels with 
only quark-antiquark interpolating fields, which may lead to unphysical results close to 
multi-hadron thresholds [12-14]. 

Here we present a lattice QCD simulation of the vector = 1 ) and scalar (0^'*') 

charmonium resonances above DD threshold, taking into account their strong decay to DD. 
The lowest vector resonance above open charm threshold, the ■0(3770) is well established 
in experiment [15], and we extract its width by simulating DD scattering in p-wave. In 
contrast to that, the experimental and theoretical status of scalar charmonia is puzzling: 
the only well-established state is the ground state Xcoi^P)^ while there is no commonly 
accepted candidate for its first excitation Xcoi'^P)- We present a study of DD scattering 
in s-wave, aiming to address this open problem. We also consider possible effects of the 
DD threshold on the vector '0(25') and scalar XcO charmonia, which lie below threshold. 

To study the hadrons present in these two channels our analysis makes a number of 
simplifying assumptions based on phenomenology, model calculations and experimental 
data: 

• We only include interpolating fields of a quark-antiquark and meson-meson type. 
For the meson meson interpolators we restrict our study to DD interpolators for the 
'I'(3770), and DD & J/'ko; interpolators for x'co- 

• We assume that elastic decay into DD is a good approximation for extracting the mass 
and width of the states. For the 'k(3770), where Hr®^P[0(377O) —>■ DD] = 93±9%, we 
neglect all further possible two- and more hadron channels, in particular decay into 
light hadrons through charm annihilation diagrams, J/Ttf-tt, J/'kr/ and ijcUi. Possible 
influence of thresholds XcO^, Xci^ and DD*, which lie above 0(3770), is also omitted. 
For the study of the scalar channel we assume elastic scattering of DD and neglect 
all other open channels such as the rjcTj or Xci^l- Note that r] is a flavour singlet on 
our Nf = 2 ensemble, and is therefore relatively heavy. 

• In the scattering analysis of the discrete energy levels we model the scattering ampli¬ 
tude using various model forms. For the 'k(3770) our model assumptions are based 
on phenomenology (experiment branching fractions and upper limits), while they are 
based on potential model expectations in the case of the Xco- 

2 Open questions for charmonia of interest 
2.1 Vector charmonia 

The 0(3770) with M = 3773.15 ± 0.33 MeV and T = 27.2 ± 1.0 MeV is located only 
~ 45 MeV above DD threshold [15, 16]. We focus here on its dominant decay mode 
0(3770) —)• DD in p-wave with branching fraction 0.93 [15]. It is a well-established 
experimental resonance and is generally accepted to be predominantly the conventional 
IDi cc state [17-19]. There is an ongoing discrepancy between results from 
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BES-II [20] and Cleo [21] regarding the non-Z)Zl part of the branching fraction which may 
be connected to neglecting interference effects in the BES-II analysis. Significant non-DD 
decays into light hadrons can occur if there is non-negligible mixing with the 'ip{2S) [ 22 ]. 
For our analysis we neglect disconnected contributions that would cause decay into light 
hadrons and treat the decay into DD as elastic, neglecting the decays into and 

J/iprj that have tiny branching fractions [15]. Our aim is to perform a determination of 
the 'ip{3770) resonance mass and ^^(3770) —)■ DD decay width using a lattice simulation for 
DD scattering in p-wave. 

We also investigate whether the DD threshold has any effect on 'ip{2S), which is the 
first radial excitation of J/V' and is situated ~ 42 MeV below threshold. Such a possibility 
was discussed in relation to the Fermilab-MILC preliminary results [23] where a simple 
analysis of the spin-averaged 25 state appeared high with respect to experiment, although 
large systematic uncertainties related to excited state contaminations were observed. A 
more recent HPQCD study [2] finds no significant discrepancy. The mixing of the vector 
charmonia with a pair of two charmed mesons was hrst simulated in [7], where only DiD 
in s-wave was considered and the width of V’(3770) was not extracted. 

2.2 Scalar charmonia 

The only well established scalar charmonium state is the ground state Xco(l.P)) interpreted 
as the ^IPq cc and located well below the open charm threshold. A further known resonance, 
the A(3915) with T = 20±5 MeV is seen only in and 77 decay channels [15]. BaBar 

has determined its quantum numbers to be O"*" [24] which would only allow = 0'*'^. 
This spin-parity determination by BaBar assumes that a = 2"*“ resonance would be 
produced in the helicity 2 state, which might not be justihed for an exotic meson^ [5]. As 
a consequence, the PDG recently assigned A(3915) to be Xco{2P) [15], but a number of 
convincing reasons given by Guo & Meissner [26] and Olsen [27] raise serious doubts about 
this assignment: 

• The dominant decay mode of scalar charmonium above open charm threshold is 
expected to be a ”fall-apart” mode into DD that would lead to a relatively broad 
resonance. In particular the width into DD is expected to be much larger than for 
the well-established Xc2{2P) [15], which decays to DD in d-wave. Yet invariant 
mass spectra of several experiments show no evidence for A(3915) —)• DD. This also 
indicates that the DD width extracted from the present lattice simulation cannot be 
compared to A(3915). 

• The spin-splitting ?tIxc 2 ( 2 P) “ ''^XcofsP) within this assignment seems too small com¬ 
pared to r?T,^^ 2 ( 2 P) - ”^X 6 o( 2 P) ”^Xc 2 (lP) “ "^Xco(lP)- 

• The partial width for the OZI suppressed A(3915) —)■ wJ/V' seems too large [26], 
which is translated to two contradicting limits for this decay in [27]. 

^For arguments in favor of the Y(3915) as a J = 2 resonance see Ref. [25]. 
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Ensemble (1) 

Ensemble (2) 

Nl X Nt 

16^ X 32 

32^ X 64 

Nf 

2 

2+1 

a [fm] 

0.1239(13) 

0.0907(13) 

L [fm] 

1.98(2) 

2.90(4) 

rriTT [MeV] 

266(3)(3) 

156(7)(2) 

LtUtt 

2.68(3) 

2.29(10) 

Kc (val) 

0.12300 

0.12686 

#configs 

279 

196 


Table 1. The gauge configurations of ensemble (1) are from [32, 33]. Those of ensemble (2) are 
provided by the PACS-CS collaboration [34]. Nl and Nt denote the number of lattice points in 
spatial and time directions, Nf the number of dynamical flavors and a the lattice spacing. 


The intriguing Xco{‘^P) was related to the broad structures in DD invariant mass in 
the same references [26, 27]. The process 77 —)• DD from BaBar [28] and Belle [29] leads 
Guo&Meissner to^ 


[26] : m = 3837.6 ± 11.5 MeV, T = 221 ± 19 MeV, (2.1) 

while e"’“e“ —)• J/xf^ DD from Belle [30] leads Olsen to 

I Q 1 ^ 

[27] : m = 3878 ± 48 MeV , T = 347 MeV. (2.2) 

Obviously the spectrum of scalar charmonia beyond the ground state presents an open 
question. Our aim is to shed some light on this issue by simulating DD scattering in s-wave 
on the lattice, and look for possible resonances in the extracted scattering matrix. Prelim¬ 
inary results based on the same simulation and only one ensemble have been presented in 
Ref. [31]. 


3 Lattice setup and charm-quark treatment 


The simulation is performed on two lattice ensembles with the parameters listed in Table 
1. Both ensembles have rather low mT^L but this is not a serious issue for charmonia and 
DD scattering in this simulation, where pions do not enter explicitly. Further details about 
the ensembles and our implementation of charm quarks may be found in [12, 14, 32, 33] 
for ensemble (1) and in [34, 35] for ensemble (2). 

To minimize heavy-quark discretization effects at finite lattice spacing the Fermilab 
method [36, 37] is used for the charm quarks. The corresponding dispersion relation [38] 
for a meson M containing charm quarks is 


EMip) = Ml + 


a^W4 


2M2 


^P' 


2\2 


(P^) 

8M| 


+ 


(3.1) 


^Here possible feed-down from 77 —>■ D*D followed by D* —>■ D'K{'y) is ignored according to [30]. 
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meson 

mass 

Ensemble (1) 

Ensemble (2) 

D 

aMi = amo 

0.9792(11) 

0.75317(83) 

D 

aM2 

1.107(14) 

0.840(22) 

D 

aMi 

1.060(44) 

0.95(15) 

spin-aver. cc 

aMi = am 

1.52451(44) 

1.20444(15) 

DD vs. cc 

2m D — El 

0.6910(36) 

0.6568(36) 


Table 2. The parameters in the dispersion relation (3.1) for D mesons and spin-averaged charmo- 
nium -|- 3mj/^). The last line is in GeV, others in lattice units. 

where p = and q G N^. 

On both ensembles the charm quark hopping parameter Kc is tuned [14, 35] using the 
spin-averaged charmonium mass 

fh = + 3mj/^) , = 3.06859(17) GeV (3-2) 

which is the relevant reference mass for our spectra of charmonium. The Mi^ 2,4 for the spin- 
averaged charmonium were determined based on the lattice data from the lattice dispersion 
relation (3.1), setting IT 4 to zero. Then Kc was fixed by tuning the kinetic mass M 2 to 
'pjjg corresponding values for the spin-averaged Mi are given in Table 2. 

To investigate the DD scattering we need the dispersion relation E£){p) for D mesons, 
which is also given by Eq. (3.1) with parameters Mi^ 2,4 in Table 2. The common feature 
of spectra in the scalar and vector charmonium channel are two-particle states DD that 
have a discrete spectrum on the finite lattice. In the absence of interactions, D{q)D{—q) 
have energies according to (3.1) 

= 2ED(q¥) , q e iVy (3.3) 

which will be shifted due to the interaction. 

Within the Fermilab approach, the rest masses have large discretization effects but 
mass differences are expected to be close to physical [39] and can be compared to experi¬ 
ment. In order to compare the splitting with we will sometimes 

plot 

E = + fh^^P (3.4) 

and compare it with E^^p _ 

An important quantity is the position of the DD threshold with respect to our reference 
mass. The splitting 2mD — Ei for ensemble (2) is very close to the experimental value 
2m^^P — frfxp ~ 0.666 GeV, while it is a bit larger for ensemble (1) due to the heavier pion 
mass and larger discretization effects (see Table 2). 

Our charm quark treatment has been verified on ensemble (1) for low-lying charmonia, 
D meson resonances [14] and Dg mesons [35, 40], where reasonable agreement with exper¬ 
iment was found. The spectrum for Dg mesons and some other hadrons containing charm 
quarks were also determined on ensemble (2) [35, 40] with even better agreement due to 
the lower pion mass and smaller discretization effects. 
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4 Analysis details 

Interpolating fields O are used to create and annihilate the physical system with = 1 
or O'*”*", isospin 1 = 0 and total momentum zero. All quark fields in the interpolators are 
smeared according to the distillation method q = Qpoint [41, 42]. We use 

Ay = 192 eigenvectors of the lattice laplacian for ensemble (2) and Ay = 96 or 64 
for ensemble (1). The distillation method is convenient for calculating a variety of Wick 
contractions. The full distillation method [41] is employed on ensemble (1) with a smaller 
volume and details of the implementation are given in [12, 14]. The stochastic version [42] 
is used on ensemble (2) with larger volume and details of our implementation are provided 
in [35]. 

4.1 Vector channel 

DD in p-wave is the dominant two-meson contribution for A < 4 GeV, while DiD appears 
higher. Sixteen cc and two DD interpolating fields are used in the relevant irreducible 
representation : 

= cAiC , (4.1) 

^15 ~ RijkC'y^E C , Ei — Qijk^~j^k , 

OTe = RijkCitl^E^c , 

Oi^ = [c75tt(ej) u'y5c{-ei) 

- c'y5u{-ei) u'y^c{ei)] -h {w d} , 

= [n5ltu{ei) U'y5'ytc{-ei) 

- C'y5ltu{-ei) il757tc(ei)] + {u ^ d} , 

where i denotes polarization, while Qijk and Rijk = Rjik are listed in [43]. The cc inter¬ 
polators Oi^i 4 for vector channel are listed in Table X of [14]. The momentum is 
projected for each D meson separately, 

iirc(k) = t)rc(x, t) , (4.2) 

X 

so that the couple to p-wave. For ensemble (1) Ay = 64 is used for , and Ay = 96 
for the remaining interpolators. 

The irreducible representation contains = 1 states of interest, and also 
"03 states with = 3 coupling due to the broken rotational symmetry on the lattice. 
In the continuum limit, contain only 1 , while contain 1 and 3 [43], 

which will help us to identify the spin 3 admixture related to 

4.2 Scalar channel 

DD in s-wave and are the dominant two-meson states in the energy region of interest 

A < 4 GeV. Seven cc, four DD, and two J/if)uj interpolating fields are used in the relevant 
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irreducible representation Af^: 

OfLr = c^c , (4.3) 

of ^ = C75tt(0) u75c(0) + {u^ d} , 
of ^ = C757iu(0) U757tc(0) + {u ^ d} , 
of ^ = X] c75M(efc) uj5c{-ek) + {u ^ d} , 

of ^ ^ c'y^u{uk) U'y^c{-Uk) + {u d} , 

|ufc|2=2 

^ ^ C7jc(0) [U7ju(0) + {u ^ d}] , 
j 

= ^C7j7tc(0) ['d7j7tu(0) + {u ^ d}] . 

j 


OfLj are listed in Table X of [14]. The momenta are projected for each meson separately 
in and . For ensemble (1) = 64 is used for Off, and iVy = 96 for 

the remaining interpolators. 

The irreducible representation Af^ contains = O^"*" states of interest, and in 
general also states with J > 4, which appear at energies beyond our interest. 

The interpolator Of^ is not used for ensemble (1) since D{2)D{—2) appears above 
4GeV. TheO-^/^^ are not used on ensemble (2) since the results from ensemble (1) indicate 
that J/V’w is almost decoupled from the rest of the system.^ 



cc 

DD 

J/'lpLO 

cc 

-o 


CM 

1 

DD 

+2 i 

+2 -4 


J/lpUJ 



+2 —4 ... 

•- --►- -• ^ 


Figure 1. Wick contractions computed for the correlation matrix (4.4) with interpolators (4.1,4.3). 
We omit contractions where the charm quark annihilates. A red solid line represents a c quark, 
while the black dashed line represents a it or d quark. 


®When the interpolators are removed from the interpolator basis, the energies En and overlaps 

{Ok\n) for the remaining eigenstates n are practically unchanged for ensemble (1). 
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4.3 Towards the spectrum 

The correlation matrix 

Cjk{t) = {a\0,{t' + t)Ol{t'm = ^ (4.4) 

n 

contains the information on energies and the overlaps Z” = {Vl\Oj\n). We evaluate 
all Wick contractions for O ~ cc, {qc){cq), {cc){qq) (4.1,4.3) shown in Fig. 1. We omit 
Wick contractions where charm quark annihilates as in almost all previous lattice simula¬ 
tion of charmonia; these induce mixing with 1 = 0 decay channels containing only light 
quarks tt, d, s, they are Okubo-Zweig-Iizuka suppressed and present a challenge for current 
lattice simulations. It is noteworthy that these decays might be important to clarify the 
experiment puzzle with regard to non-DD hadronic decays [18, 22]. 

The energies and overlaps are extracted from the correlation matrix using the gener¬ 
alized eigenvalue method [44-47] 

{t) = {t)C{to)u^^'> {t) , (4.5) 

where A*'”'^(t) oc at large t. Correlated two or one-exponential fits to X^"'\t) are 

used and to = 2, 3. The errors-bars correspond to statistical errors obtained using single¬ 
elimination jack-knife. 


ensemble (1) ensemble (2) exp. 



(a) (b) (c) (d) (e) 


only cc only cc 

Figure 2. The energies E (see Eq. (3.4)) in the vector channel on both ensembles, together with the 
experimental masses. The circles represent = 1 states, while the diamond represents a 3 
admixture present in the irreducible representation and related to the tps. The dashed lines 
show the non-interacting energy of D{1)D{—1) (3.3), and the dotted line represents the threshold 
2mD- The D{0)D{0) state does not appear for p-wave. Interpolators used in (a,c) are given in 
Table 3, while (b,d) utilize just from the same sets. 

























n 

fit 

range 

fit 

type 

x" 

d.o.f. 

Eiata 

E [GeV] 
(3.4) 


Ens. 

(1) 




1 

3-14 

2e^ 

8.57/8 

1.54153(43) 

3.09572(34) 

2 

3-14 

26^= 

16.62/8 

1.9045(38) 

3.6738(58) 

3 

3-13 

2e= 

5.18/7 

1.9801(46) 

3.7941(71) 

4 

3-13 

2e= 

5.09/7 

2.0109(60) 

3.8433(93) 

5 

3-13 

2e= 

8.49/7 

2.1105(21) 

4.0019(32) 


Ens. 

(2) 




1 

3-29 

2e= 

3.15/23 

1.21683(16) 

3.09557(18) 

2 

3-11 

26^= 

3.44/5 

1.4862(60) 

3.682(13) 

3 

3-11 

2e‘= 

4.36/5 

1.531(11) 

3.779(24) 

4 

3-11 

2e‘= 

4.92/5 

1.5611(78) 

3.845(17) 

5 

3-11 

2e^ 

4.78/5 

1.5661(75) 

3.856(16) 


{apf 

(ap)® cot(5) 

(ap)^ cot(5) 

0°1 

/ 

/ 

/ 

/ 

-0.0588(47) 

0.0137(18) 

0.00717(95) 

-109.1(6.4)i 

0.02413(57) 

-0.00599(34) 

-0.00303(17) 

148.0(7.7) 

/ 

/ 

/ 

/ 

0.1755(33) 

-0.144(23) 

-0.068(11) 

153.0(4.4) 

/ 

/ 

/ 

/ 

-0.0169(50) 

0.0021(10) 

0.00143(68) 

120(25)i 

0.0207(93) 

0.00056(255) 

0.00037(167) 

79(40) 

/ 

/ 

/ 

/ 

0.0509(65) 

-0.0054(76) 

-0.0034(49) 

115(35) 


Table 3. Discrete lattice spectrum from charmonium in the irreducible representation which 
contains = 1 , 3 and higher J states. The p and 5 correspond to DD scattering in 

p-wave. Subset from the interpolators in Eq. (4.1) is used for ensemble (1) 

and 23.15) ^17^ fo'' ensemble ( 2 ). to = 2 is used for all data points. 


ground st. 

Ensemble (1) Ensemble (2) 

Eiata 

2.0124(38) 1.559(51) 


Table 4. The energy of ipz with = 3 from the ground state in the A 2 irreducible 
representation of the Oh point group. 

5 Results for the vector channel 

5.1 Discrete spectrum 

The energy levels in the vector channel are shown in Fig. 2a and 2c together with the 
experimental masses. The full set of operators gave noisier signals than suitable subsets, 
and the chosen subsets are listed in Table 3. The circles denote the energy levels that are 
related to = 1 states J/V', '0(25'), 0^(3770), D{1)D{—1) (from bottom to top), 
while D{0)D{0) does not appear for p-wave. The diamond indicates a level related to the 
= 3 state 0^3, that is present in representation due to the broken rotational 
symmetry on the lattice. 

The highest state {n = 5) has largest overlap with and disappears when these 

interpolators are excluded from the basis, as shown in Fig. 2b and 2d. Each energy level 
in addition to D{1)D{—1) indicates the presence of a bound state or a resonance. Good 
resemblance with the experimental spectrum is indeed confirmed in Fig. 2. The J/ij: is 
significantly below threshold and no effect from threshold is expected. The 0(25) is situated 
~ 42 MeV below threshold in experiment, and the corresponding finite volume energy on 
the lattice does not depend (within uncertainties) on whether DD interpolators are used or 
not (see Fig. 2). The appearance of levels n = 3 and 4 is related to the 0(3770) resonance 
and to the spin 3 admixture and the corresponding 03 resonance. Level n = 4 is related 
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Figure 3. The overlaps Z” = {U\Oj\n) for the vector channel show the matrix elements of 
interpolators Oj between the vacuum (n| and the eigenstate \n) on the lattice. We present the 
overlap ratios fmaxmZ'^ on ensemble (1) (top) and on ensemble (2) (bottom). The denominator 
is the maximal \Z‘^\ at given operator number j. These ratios are independent on the normalization 
of the interpolators Oj. Levels n = 1,.., 5 are ordered from lowest to highest En in Figs. 2a and 2c 
for both ensembles, respectively. The order of interpolators j on the abscissa (listed in caption of 
Fig. 3) is the same as in the list (4.1). 


to ^jJ 3 due to smaller overlaps = 4). This is based on the fact that Oi(li 4 couple 

in the continuum limit only to 1 (which is responsible for small (Oi(Li 4 |V' 3 ) s-t finite 
a), while Ofl couple to 1 and 3 . Further support is given by the near-degeneracy 

with the energies from the irreducible representation ^4^” where a 3 state comes as the 
ground state (see Table 4). For the ^^(3770) the avoided-level crossing scenario suggests 
in the energy region m ± F, which is reasonably satisfied by comparing to experiment. In 
order to really determine the resonance mass and width for 'ip{3770) one needs to consider 
the phase shifts for DD scattering in p-wave. 


5.2 DD scattering in p-wave 

We assume that DD scattering in p-wave near the resonance V’(3770) is elastic, which is a 
good approximation since Br[ip{3770) —)• DD] = 93 ± 9%, while the remaining part goes 
mainly to light hadrons and charmonium states (i.e. J/ipnir, J/iljr ]..In the elastic case, 
the scattering phase shift 5 is given by Liischer’s relation [48, 49] 


p cot 5 (p) 


2Zoo(l;(g§)^) 

Ly/ir 


(5.1) 


which applies for the total momentum zero employed in our case. The momentum p oi D 
mesons is extracted from the measured energy levels = 2E£){p) using the dispersion 
relation (3.1). The resulting momenta and phase shifts for all eigenstates except for the 
spin 3 admixture and for the finite volume state related to J/'ij: are collected in Table 3. 
The large absolute value of p^ cot 5 corresponds to feeble scattering, while small p^ cot 5 is 
related to significant scattering. 
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We fit our data in two ways: 


ensemble (1) 


ensemble (2) 




I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I I 

-0.20 -0.10 0.00 0.10 0.20 0.30 0.40 0.50 

p^[GeV^] 


Figure 4. cot S /i/s versus for DD scattering in p-wave in the region of the bound state 

and the '0(3770) resonance. The p denotes the momentum of D meson. We show the Breit-Wigner 
fit (i) and the extended fit (ii), which aims to capture also the behavior around ijj(2S). 



Ensemble (1) 
fit (i) fit (ii) 

Ensemble (2) 
fit (i) fit (ii) 

exp 

D+D-/D^D^ 

0(3770) 
PR [GeV] 
mu [GeV] 
g (no unit) 

0.208(31)(3) 0.159(35)(2) 
3.784(7)(10) 3.774(6)(10) 
13.2(1.2) 19.7(1.4) 

0.334(155)(5) 0.343(184)(5) 
3.786(56)(10) 3.789(68)(10) 
24(19) 28(21) 

0.26/0.29 

3.77315(33) 

18.7(1.4) 

0(25) 
\pb\ [GeV] 
niB [GeV] 

0.380(17)(6) 

3.676(6)(9) 

0.280(43) (4) 
3.682(13)(9) 

0.31/0.28 

3.6861090^4 


Table 5. Parameters of the resonance '0(3770) and bound state '0(25') from fits (i) (5.4) and 
(ii) (5.6). The 0(3770) —>■ DD width T = g^p^ / {Qtts) is parametrized in terms of the coupling 
g and compared the value of the coupling derived from experiment [15]. The pr denotes D- 
meson momenta at the peak of the resonance and \pb\ the binding momentum. The first errors 
are statistical and the second errors (where present) are from the scale setting uncertainty. The 
experimental data and errors are based on PDG values. Errors on experimental Pr/b are suppressed 
as they are very small. 

(i) A resonance 0(3770): The scattering matrix in the vicinity of a resonance has a 
Breit-Wigner form 


The width 


_ \/sr(s) _ 1 

— s — iyfs T{s) coi 6i{s) — i 


(5.2) 


(5.3) 
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is parametrized in terms of the phase space for p-wave decay and the "0(3770) —)• DD 
coupling g. It is expected that the leading dependence of F on rriu/d is captured by phase 
space. Equations (5.2, 5.3) lead to cot (5i(s)/y^ = (67r/p^)(m^ — s) and then expressing 
s = 2(m|, +p^)^/^ and niR = 2(m^ + to 


p^ cot (5(s) 


Gvr 




(5.4) 


where pR is the D meson momentum at the resonance peak. The values of g and pR follow 
from the linear fit (5.4) through the energy levels n = 3,5 in the vicinity of the resonance, 
where the Breit-Wigner form applies (level 4 is omitted since it is attributed to 03 as 
discussed above). The fit is shown in Fig. 4, while the resulting resonance parameters are 
given in Table 5. The resonance mass niR corresponding to the pR on the lattice is given 
by inserting the Fermilab dispersion relation (3.1) in (3.4) 


rriR/B = 2 Ed{,Pr/b) - 


(5.5) 


and will be used for resonances or bound states throughout this work. 


(ii) A resonance 0(3770) and a bound state 0(2S): In addition to the Breit-Wigner 
form (5.4), which is linear in p^, we make use also of the square form in p^ 


V 


cot (5i (s) = A -|- Bp^ -|- Cp^ 


(5.6) 


which in general has a longer range of applicability. It aims to capture also the DD 
scattering in the vicinity of 0(25'): there the (imaginary) phase shift in Table 3 nearly sat¬ 
isfies the condition for the bound state cot <5 ~ z on the physical Riemann sheet pb = 
i\pB\, leading to p^cotJ ~ \pb\^- The fit (5.6) through levels n = 2,3,5 in Fig. 4 
assumes that the 0(25) state still affects the DD scattering. It renders (A, B, C) ~ 
(0.0046(19) GeV0-0.168(27),-0.52(13)/GeV2), (0.0069(88) GeV0-0.023(80),-0.30(68)/GeV^) 
on ensemble (1) and ensemble (2) respectively. The zero ofp^cot(ii/-0s, and the derivative 
at this zero, lead to the parameters of 0(3770) resonance in Table 5. This model also leads 
to a bound state 0(25) at pB = i\pB\ where the scattering amplitude T (5.2) has a pole 
and cot 6{pb) = i- The bound state mass tub in Table 5 is indeed close to experimentally 
measured 0(25). 

On ensemble (2) the results both from fit (i) and fit (ii) are compatible with the 
experimental data'^ within large statistical uncertainties (see Table 5). Note that the higher- 
lying 0(4040) resonance does not influence the results (for this ensemble), since it lies 
significantly higher than the relevant energy levels. 

“^Since we work in the isospin-symmetric limit we measure the sum of the neutral and charged decay 
modes; therefore we compare to the experimental value + '?d+_d- obtained from r[i/)(3770) — >■ 

DoDo] = g^j^iQBoP^/{tms) and r0(377O) —>■ D'^D~] — [15]. Notice also that averaging 

the results from recent experiment resonance mass determinations for the ip{3770) leads to a value of 
— 3778.1(1.2), much larger than the fit by the PDG (which relies on an experiment neglecting 
interference with non-resonant background) and consistent with the most recent results in [16]. 
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On ensemble (1) the results for 'ip{3770) from fit (i) give a smaller resonance momentum 
Pr than in experiment, which we attribute to the unphysical threshold on ensemble (1) at 
ttItt — 266 MeV and the finite lattice spacing. The resonance mass rriR calculated as in 
Eq. 3.4 compares favorably. The coupling constant from fit (i) is to small compared to 
experiment which is likely related to to the closeness of the '0(4040) resonance neglected 
in the analysis. The assumption that the resonance ■0(4040) does not affect the energy 
level related to D{1)D{—1) is probably not justified on ensemble (1), where energy level 
lies higher (and closer to 0(4040)) than on ensemble (2). Roughly estimating the effect by 
comparing the one-resonance and two-resonance scenarios, estimating g andpR for 0(3770) 
and 0(4040) from available experimental data [50], the coupling we observe is consistent 
with this interpretation ^ . Given the possibly large influence from the 0(4040) we can not 
conclude that ht (ii) is better than ht (i) on this ensemble. 

The resulting 1 spectrum is summarized and compared to experiment in Fig. 5. 


ensemble (1) ensemble (2) exp. 



fit (i) fit (ii) fit (i) fit (ii) 


Figure 5. The comparison of the final 1 spectrum to the experiment. The magenta diamond 
denotes 0(3770) resonance mass from the Breit-Wigner fit (i) or extended fit (ii), given in Eqs. (5.4) 
and (5.6), respectively. The magenta triangle denotes 0(25') obtained as a pole in DD channel. 
The blue triangles denote masses of J/0 and 0(25) extracted as energy levels in the finite box. 
The statistical and scale setting errors have been summed in quadrature. 


®The maximal effect of 0(4040) is estimated by assuming that 0(4040) width is saturated by DD (instead 
of DD, DD*, D*D* and other modes). 
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Figure 6. The energies E (see Eq. 3.4) in the scalar channel on both ensembles. The only 
well established experimental state Xco(l^’) is shown by the (magenta) circle. Triangles show 
three intriguing candidates for Xco(2T’), that are not universally accepted: 2f(3915) and the broad 
resonances (2.1,2.2) suggested in [26, 27]. The dashed lines shown energies of non-interacting 
D{q)D{—q) with q = 0,1,2 (3.3), while dot-dashed line represents rrij^^ + Interpolators used 
in (a,d) are given in Table 6, (b) uses in addition, while (c,e) are based only on 5. 


6 Results for the scalar channel 
6.1 Discrete spectrum 

The energy levels in the scalar channel are shown in Figs. 6. The only experimentally well 
established state is Xcoi^P)- The triangles represent the intriguing experimental candidates 
for Xco(2T’), none of which is commonly accepted (see Section 2). 

The spectrum from a lattice simulation consists both of energy levels that have large 
overlap with qq operators as well as energy levels with dominant overlap to DD operators. 
The latter appear near their non-interacting energies of Eq. (3.3), which are denoted 
by dashed lines in Figs. 6a, d. On ensemble (1) levels n = 2,4 appear near the non¬ 
interacting D{0)D{0) and D{1)D{—1) (cf. Fig. 6a). Levels n = 2, 3,4 on ensemble (2) have 
dominant overlap to DD scattering operators and are close to non-interacting D(0)Z)(0), 
D{1)D{—1) and D{2)D{—2) energies (cf. Fig. 6d). 

In the elastic case each energy level in addition to the number of expected D{q)D{—q) 
scattering levels is related to the presence of a bound state or a resonance. There are two 
such states, that cannot be attributed to D{q)D{—q) for both ensembles. The ground state 
is related to Xcoi^P) and is close to its experimental mass. The second of these two levels 
appears above threshold and corresponds to n = 3 for ensemble (1) and n = 5 for ensemble 
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n 

fit 

fit 

d.o.f. 


E [GeV] 

{apf 

(ap) cot ((5) 

(ap) cot(5) 



range 

type 



(3.4) 






Ens. 

(1) 








1 

6-15 

2e= 

9.50/6 

1.7468(19) 

3.4226(27) 

-0.2226(49) 

-0.4716(52) 

-0.2700(31) 

-240.9(2.8)1 

2 

6-15 

2e= 

5.40/6 

1.9494(33) 

3.7453(52) 

-0.0099(32) 

0.11(11) 

0.058(55) 

4(190) + 84(306)1 

3 

3-12 

2e= 

2.73/6 

2.0625(81) 

3.925(13) 

0.1185(98) 

0.39(17) 

0.191(79) 

41(11) 

4 

3-12 

2e= 

6.39/6 

2.1190(31) 

4.0154(49) 

0.1857(48) 

-0.50(11) 

-0.236(52) 

139.2(6.5) 


Ens. 

(2) 








1 

4-15 

26^= 

7.84/8 

1.3672(38) 

3.4227(83) 

-0.1136(60) 

-0.3371(89) 

-0.2466(68) 

-344.9(8.8)1 

2 

4-15 

26^= 

3.74/8 

1.5018(74) 

3.715(16) 

-0.0038(61) 

-0.004(185) 

-0.003(123) 

-4(201)1 

3 

3-12 

26^= 

3.50/6 

1.5497(71) 

3.820(16) 

0.0367(61) 

1.2(7.5) 

0.8(4.8) 

8.7(29.2) 

4 

6-12 

le= 

0.89/5 

1.5934(92) 

3.915(20) 

0.0745(80) 

1.7(10.8) 

1.1(6.8) 

9.1(28.3) 

5 

6-12 

le= 

1.68/5 

1.6148(85) 

3.961(19) 

0.0932(75) 

-0.21(21) 

-0.13(13) 

124(28) 


Table 6. Discrete lattice spectrum for the scalar channel. The p and 5 correspond to DD scat¬ 
tering in s-wave. Subset 0^3 5 ,from interpolators in Eq. (4.1) is used for ensemble (1) and 
0^3 5 , for ensemble ( 2 ). to = 2 is used for all data points. 


(2), as shown in Figs. 6a and 6d. The avoided level crossing scenario suggests that an 
additional level appears somewhere in the range ~ m ± F, which suggests the existence 
of a resonance roughly at 

m ~ 3.9 — 4.0 GeV (naive estimate from En) . (6-1) 

This is close to the first excitation obtained using just 0“ interpolators in Figs. 6c and 6e. 
Such a basis gives a rough estimate of resonance masses but is not well suited to capture 
two-particle states or resonances and bound states close to threshold [12-14]. 

The spectrum including J/'tpuj interpolating fields is shown in Fig. 6b for ensemble (1). 
An energy level related to J/'ip{0) w(0) appears at roughly + while the energies of 
all the other levels remain unaffected with respect to Fig. 6a. We have verified also that 
the overlaps for the remaining levels are not affected if are in the basis or not. This 

indicates that the J/'tpuj channel is decoupled from DD channel to a good approximation. 

6.2 DD scattering in s-wave 

We now study DD scattering assuming that J/'tjjuj channel is decoupled as argued in 
the previous paragraph. We did not include r]cf] interpolators in the correlation matrix 
assuming that they can be neglected. The energy shifts of the extracted with respect 
to E/j/o give the size of the s-wave scattering phase shift 5 according to (5.1). On ensemble 
(1) we observe statistically significant energy shifts with respect to the dashed lines in Fig. 
6 a. The energies yield T)-meson momenta p via El//*' = 2E]p{p), and the corresponding 
phase shifts 6{p) via Eq. (5.1). These are provided for all levels in Table 6 and plotted in 
Figs. 7 and 8. 

The uncertainties on the energies E'„= 2 , 3,4 are rather large for ensemble (2) and they are 
within errors compatible with non-interacting energies E/j/^ (3.3). This implies that we are 
not able to reliably determine the energy shifts, and the resulting errors on the scattering 
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Figure 7. We show p cot 6 and p cot 6/y/s versus p^ for DD scattering in s-wave, where p denotes 
the momentum of the D meson. The circles with (sizable) errors denote the lattice data, while the 
solid lines show pcoi5 = 2ZqqI{L^Jtt) according to Liischer’s relation (5.1). When the momentum 
is compatible with the non-interacting momentum p = 27rq/L (q G N^), one has 5 = 0 and 
I cot 5| = oo, which is responsible for the huge errors on pcot 6 on ensemble (2). 


matrix will be large, as illustrated in Fig. 7. If ~ within errors, this implies 

(5 ~ 0 modulo vr and cot 5 ~ ±oo within errors. The extracted pcotd from n = 2,3,4 
have large errors, which allow almost all pcot6 expect for small |pcot5|. For n = 2,3,4 
we plot central values with the ranges [f{p^ — o'p 2 ),/(p^ -|- 0 -^ 2 )] where / = pcot5 

or pcot5/y/s. The error for all other levels (on both ensembles and both channels) is the 
usual jack-knife. 

The resulting pcot6/^/s for ensemble (1) has a puzzling behavior and we are going 
to confront it with various hypothesis, collected in Fig. 8. The errors of pcot6/y/s on 
ensemble (2) are large and do not allow reliable fits. We will still compare data from 
ensemble (2) with fits that are based on ensemble (1), and plot them as function of p^ on 
the same figure. 

(i) A narrow resonance: In the vicinity of a Breit-Wigner resonance (5.2) one expects 


pcot 5{s) 


X{pr-p^): ^{s) = g^- 

9 ^ 


( 6 . 2 ) 


and the zero gives the position of the resonance. The upper three points in Fig. 8a however 
do not fall onto one line, so our results cannot be reconciled with a single Breit-Wigner 
resonance in the region between 2m d and 4 GeV. The highest two points support the 
existence of a narrow resonance between them and a linear fit (6.2) over two levels shown 
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Figure 8. The p cot S / y/s versus for DD scattering in s-wave, where p denotes the momentum 
of the D meson. The lattice data (blue and red circles) is confronted with p cot S/,/s based on 
various hypothesis (dashed lines) described in Section 6 of the main text. The thin dot-dashed lines 
in the plots at the bottom denote pcotS = [L/H') (5.1). The left-hand column shows the 

results for ensemble (1), the middle column for ensemble (2) and the left-hand column is an overlay 
of both for comparison. 

in Fig. 8(i,a) renders (5.5) 


mn = 3.966(20) GeV, g = 1.26(18) GeV, (6.3) 

= 0.614(33) GeV, = 62(17) MeV . 

The scattering data from ensemble (1) therefore suggest the existence of a yet unobserved 
scalar state called x'co- Note that its mass is within the range of the naive estimate (6.1). 
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It has a width in our simulation, while the corresponding width in experiment 

would be modified due to a different phase space via p = [m^/4 — leading to 


^7ed^ct = 67(18) MeV . 


(6.4) 


We have assumed that g and rriR do not depend on the pion mass here. It is unlikely that 
this state corresponds to X(3915) since the DD decay channel was not observed for this 
state. A narrow resonance is roughly consistent also with the result from ensemble (2) 
within huge errors (see Fig. 8(i,b)), however there must be some additional interaction 
between D and D near the threshold according to ensemble (1). 


(ii) A narrow resonance and a bonnd state Xco(lP): Our next hypothesis assumes 
that Xco(lP) represents a pole in DD scattering on the first Riemann sheet, leading to 
pcoi5 ~ i\pB\i = —\pb\ at the position of the bound state. The negative value of pcotS 
below threshold might be a possible reason why pcot6 at threshold is smaller than ex¬ 
pected based on narrow resonance (6.2,6.3). In this case the value of pcot 6 at threshold is 
influenced by the resonance and a bound state. To investigate this situation, we attempted 
several hts over all four levels on ensemble (I). A form that is motivated by our data 


pcot 5{s) 


= A + B p"^ + 


C 


(6.5) 


^/s ' ^ ' p^ — D 

is presented in Fig. 8(ii,a), where A = 0.13(15), B = 0.66(18)/GeV^, C = 0.028(63) GeV^ 
and D = 0.513(77) GeV^ are obtained from the fit. This hypothesis supports a bound 
state at ps = i\pB\ which corresponds to a pole in T (5.2) or equivalently cot 5{pb) = i, 

i.e. 

\pb\ = 0.7517(83) GeV tub = 3.4224(27) GeV . (6.6) 


The bound state is attributed to Xco(lP) and its mass is very close to the one obtained 
from the ground state energy. The hypothesis also supports a narrow resonance at pB = 
0.668(35) GeV where function (6.5) crosses zero, and 


ruR = 4.002(24) GeV, g = 0.85(65) GeV, (6.7) 

= 30(45) MeV = 32(48) MeV . 

This is roughly consistent with the x(.q in (6.3). This hypothesis based on ensemble (1) 
is consistent also with the result from ensemble (2) within huge errors in Fig. 8(ii,b). An 
interesting feature of this hypothesis is the large p cot 6 or equivalently small cross-section 
at p^ ~ D, which corresponds to ^/s ~ 4.0 GeV. This feature seems to be present also in 
the experimental data from Belle [27] where a dip seems to appear at similar invariant mass. 


(iii) A broad resonance: The broad resonances (2.1,2. 2) proposed by Meissner&:Guo [26] 
or Olsen [27] are compared with our lattice data in Fig. 8(iii). This shows a Breit-Wigner 
shape (6.2) with pR and g extracted from the experimental data (2.1,2. 2). Although they 
are roughly compatible with our scattering results near threshold, they cannot be reconciled 
with it in the region above threshold where our data indicates either a much narrower 
resonance or a more complicated situation not covered by our assumptions. 
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(iv) Two resonances: Since neither one narrow or one broad resonance describe our 
scattering data near and above threshold, we next try an hypothesis with two elastic 
resonances 


pcot (5(s) 

^ 7 ^ 


9a ^ 9l 


HPrb-p'^) 


-1 


( 6 . 8 ) 


With this parametrisation there are two resonance poles in the scattering amplitude, sep¬ 
arated by a zero. Figure 8(iv,a) shows an example with qa = 2.1 GeV, = 0.23 GeV, 
qb = 1.0 GeV and pRg = 0.57 GeV that is consistent with the upper three scatter¬ 
ing points for ensemble (1)®. This hypothesis however predicts another energy level near 
p^ ~ 0.1 GeV^ where the model (6.8) crosses with the Liischer curve. Another energy 
level is expected in the two-resonance scenario also according to naive reasoning that each 
resonance or bound state leads to a level in addition to DD. Such an additional energy 
level at ~ 0.1 GeV^ is not observed in ensemble (1) indicating that this hypothesis 
is not supported by our data. An analogous conclusion is reached when confronting this 
hypothesis with the data from ensemble (2); the hypothesis predicts five energy levels in 
the region = [—0.1,0.5] GeV and we observe four levels only. 


7 Conclusions and outlook 

We performed a lattice QGD simulation of DD scattering in s-wave and p-wave to study 
vector and scalar charmonium resonances on two rather different ensembles. This is an 
exploratory simulation and the first step towards determining the strong decay width of 
charmonium resonances above open charm threshold. Ensemble (1) has Nf = 2 and = 
266 MeV, while ensemble (2) has Nf = 2 + 1 and = 156 MeV. Several cc and DD 
interpolating fields were used in both channels, where the (stochastic) distillation method 
was used to evaluate the Wick contractions. Our analysis relies on the assumption that 
looking at elastic scattering in a single channel (DD) is a good approximation. 

In the vector channel, the well known ■0(3770) resonance is present just above DD 
threshold with ilr®^^'[0(377O) —)• DD] = 93 ± 9%. We assume that the DD scattering 
is elastic in this energy region and determine the phase shift for DD scattering in p-wave 
using the Liischer formalism. The Breit-Wigner fit is performed in vicinity of the 0(3770) to 
obtain its resonance mass at 3.784(7) (10) GeV and 3.786(56) (10) GeV for ensembles (1) and 
(2), respectively. Our determination of its decay width might be affected by the 4^(4040) 
on ensemble (1). Ensemble (2) does not suffer from this issue, and the determination 
of the resonance parameters is more reliable, but its statistical accuracy is poor. The 
resulting spectrum in the vector channel, including also J/0 and 0(25"), is compared to 
experiment in Figure 5. This work presents a step towards a determination of the 0(3770) 
resonance parameters from lattice QGD. Improvement of the results for this resonance in 
future lattice studies will need consideration of multiple volumes and momentum frames, 
further scattering channels and higher statistics. 

® These values are not obtained from a fit, bnt present one example of four parameters, where (6.8) is 
consistent with upper three scattering points. 
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In the scalar channel, only the ground state Xcoi^P) is understood and there is no 
commonly accepted candidate for its first excitation Xco{‘^P)- Guo & Meissner [26] as well 
as Olsen [27] argued that the higher lying X(3915) can probably not be identihed with 
the Xco(2-P). They suggest that a broad structure observed in the DD invariant mass 
represents Xco(2-P). This posed a particular motivation to extract the phase shift for DD 
scattering in s-wave in the present work. The resulting scattering data on the ensemble 
with ttItt = 156 MeV is unfortunately noisy. The simulation at m-,^ = 266 MeV renders the 
scattering phase shift only at a few values of the DD invariant mass, which also does not 
allow a clear answer to the puzzles in this channel. We obtain the Xco(lT’) and our data 
provides an indication for a yet-unobserved narrow resonance slightly below 4 GeV with 
r[Xco DD] below 100 MeV. A scenario with this narrow resonance and a pole in the 
DD scattering matrix at Xco(lT’) agrees with the energy-dependence of our phase shift. 
We discussed three other scenarios; just one narrow resonance, just one broad resonance 
(proposed in Guo & Meissner [26] and Olsen [27]), or one narrow and one broad resonance. 
None of these scenarios agree with our current data in the whole energy region probed, 
however we can not currently exclude these possibilities. For the scalar channel this leaves 
us with a situation where puzzles remain, both from theory and experiment. To clarify the 
situation, further experimental and lattice QGD efforts are required to map out the s-wave 
DD scattering in more detail. Including further coupled channels in the future would be 
useful to relax the model assumptions made in our current study. 
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